
model{

#== pop 1 ==#

for(i in 1:n1){
for (k in 1:K){
s1_1[i,k]<- se[k]^x1[i,k]*((1-se[k])^(1-x1[i,k]))
s2_1[i,k]<- sp[k]^(1-x1[i,k])*((1-sp[k])^x1[i,k])
}
prob1[i]=pi1*(prod(s1_1[i,1:K])) + (1-pi1)*(prod(s2_1[i, 1:K]))
z1[i] ~ dpois( - log(prob1[i]))
}

#== pop 2 ==#

for(i in 1:n2){
for (k in 1:K){
s1_2[i,k]<- se[k]^x2[i,k]*((1-se[k])^(1-x2[i,k]))
s2_2[i,k]<- sp[k]^(1-x2[i,k])*((1-sp[k])^x2[i,k])
}
prob2[i]=pi2*(prod(s1_2[i,1:K])) + (1-pi2)*(prod(s2_2[i, 1:K]))
z2[i] ~ dpois( - log(prob2[i]))
}

for (k in 1:K){
se[k] ~ dbeta( omega1*(kappa1 -2)+1, (1-omega1)*(kappa1-2) +1)
sp[k] ~ dbeta( omega2*(kappa2 -2)+1, (1-omega2)*(kappa2-2) +1)
}


omega1 ~ dbeta(1,1)T(0.5,)
omega2 ~ dbeta(1,1)T(0.5,)
kappa1 = kappaMinusTwo1 +2
kappaMinusTwo1~ dgamma(0.01,0.01)
kappa2 = kappaMinusTwo2+2
kappaMinusTwo2 ~ dgamma(0.01,0.01)
pi1 ~ dbeta(1.8,2)
pi2 ~ dbeta(1.8,2)
}
